library(copula)
library(lattice)
# Exponential distribution for both marginals
lambda <- 1
# Parameters for the marginals
# You can adjust these parameters as needed
params1 <- lambda # Rate for the first exponential distribution
params2 <- 2*lambda # Rate for the second exponential distributionModelos - Punto Extra
Intrucciones
Considera dos variables \(X\) y \(Y\) con distribución exponencial donde una tiene valor esperado que sea el doble del valor esperado de la otra (tú decides el parámetro).
- Escoge al menos 5 tipos de cópulas para modelar la dependencia entre ellas, para cada cópula escoge parámetros que reflejen una dependencia negativa y una dependencia positiva (en total tendrás 10 modelos). Para cada uno de los 10 modelos:
- Muestra la función de distribución conjunta de \(X\) y \(Y\) (la fórmula).
- Grafica la función de densidad conjunta de \(X\) y \(Y\).
- Varianza de \(X\), \(Y\) y \(X + Y\).
- Desviación estándar de \(X\), \(Y\) y \(X + Y\).
- Cuantil al 95% de \(X\), \(Y\) y \(X + Y\).
- Valor esperado condicional de la cola al 95% de \(X\), \(Y\) y \(X + Y\).
- Repite el inciso anterior; pero ahora considera que \(X\) y \(Y\) tienen distribución pareto donde una tiene valor esperado que sea el doble del valor esperado de la otra (tú decides los parámetros).
Solución
Caso Exponencial
Tenemos las siguientes variables aleatorias:
\[X\sim \textbf{Exp}(\lambda)\]
\[Y\sim \textbf{Exp}(2\lambda)\]
donde \(F_{x}(x)=1-e^{-\lambda x}\) y \(F_{y}(y)=1-e^{-2\lambda y}\)
Cópula Gaussiana
\[C(u,v)=\Phi_{2}(\Phi^{-1}(u),\Phi^{-1}(v), \rho)\]
donde \(\Phi_{2}\) es la función de distribución conjunta de una distribución normal bivariada con media cero y varianza uno y covarianza \(\rho\) (esta determinará la relación negativa o positiva de las variables); y \(\Phi\) es la función de distribución de una distribución normal estándar.
Función de distribución conjunta
\[H(x, y) = \Phi_{2}\left(\Phi^{-1}(1 - e^{-\lambda x}), \Phi^{-1}(1 - e^{-2\lambda y}); \rho \right)\]
Función de densidad conjunta
library(copula)
rho =.9
# Definir cópula gaussiana
myCopulaP <- normalCopula(rho, dim = 2)
myCopulaN <- normalCopula(-rho, dim = 2)
# Número de muestras
n <- 10000
# Generar muestras de la cópula
set.seed(123)
samples <- rCopula(n, myCopulaP) # caso positivo
samples2 <- rCopula(n, myCopulaN) # caso negativo
# Convertir las muestras a distribuciones exponenciales
## caso positivo
samples[,1] <- qexp(samples[,1], rate = params1)
samples[,2] <- qexp(samples[,2], rate = params2)
## caso negativo
samples2[,1] <- qexp(samples2[,1], rate = params1)
samples2[,2] <- qexp(samples2[,2], rate = params2)
library(ks)
# Estimación de la densidad
dens <- kde(as.matrix(samples)) # caso positivo
dens2 <- kde(as.matrix(samples2)) # caso negativo
library(plotly)
# Convertir la estimación de densidad a un formato adecuado
## caso positivo
x <- dens$eval.points[[1]]
y <- dens$eval.points[[2]]
z <- matrix(dens$estimate, nrow = length(x), ncol = length(y))
## caso negativo
x2 <- dens2$eval.points[[1]]
y2 <- dens2$eval.points[[2]]
z2 <- matrix(dens2$estimate, nrow = length(x2), ncol = length(y2))Caso Positivo
# Crear gráfico 3D
## caso positivo
plot_ly(x = ~x, y = ~y, z = ~z, type = "surface")Caso Negativo
plot_ly(x= ~x2, y = ~y2, z = ~z2, type = "surface")Varianza
# Calcula las estadísticas para las muestras generadas
# Para la correlación negativa (o positiva, según corresponda)
## Caso positivo
X_Positivo <- samples[,1]
Y_Positivo <- samples[,2]
Sum_XY_Positivo <- X_Positivo + Y_Positivo
## Caso negativo
X_Negative <- samples2[,1]
Y_Negative <- samples2[,2]
Sum_XY_Negative <- X_Negative + Y_Negative
# Varianza
var_X <- 1
var_Y <- 1/4
## Caso Positivo
var_Sum_XY_Positivo <- var(Sum_XY_Positivo)
## Caso Negativo
var_Sum_XY_Negative <- var(Sum_XY_Negative)
# Desviación estándar
sd_X <- sqrt(var_X)
sd_Y <- sqrt(var_Y)
## Caso Positivo
sd_Sum_XY_Positivo <- sd(Sum_XY_Positivo)
## Caso Negativo
sd_Sum_XY_Negative <- sd(Sum_XY_Negative)
# Cuantil al 95%
quantile_X_95 <- qexp(0.95, rate = params1)
quantile_Y_95 <- qexp(0.95, rate = params2)
## Caso positivo
quantile_Sum_XY_95_Positivo <- quantile(Sum_XY_Positivo, 0.95)
## Caso Negativo
quantile_Sum_XY_95_Negative <- quantile(Sum_XY_Negative, 0.95)
# Valor esperado condicional de la cola al 95%
## Caso positivo
#exp_X_Tail_95 <-
#exp_Y_Tail_95 <-
mean_Sum_XY_Tail_95_Positivo <- mean(Sum_XY_Positivo[Sum_XY_Positivo > quantile_Sum_XY_95_Positivo])
## Caso Negativo
mean_Sum_XY_Tail_95_Negative <- mean(Sum_XY_Negative[Sum_XY_Negative > quantile_Sum_XY_95_Negative])
# Imprimir los resultados
cat("Correlación Negativa:\n")Correlación Negativa:
cat("Varianzas: X =", var_X, ", Y =", var_Y, ", X + Y =", var_Sum_XY_Negative, "\n")Varianzas: X = 1 , Y = 0.25 , X + Y = 0.634196
cat("Desviaciones Estándar: X =", sd_X, ", Y =", sd_Y, ", X + Y =", sd_Sum_XY_Negative, "\n")Desviaciones Estándar: X = 1 , Y = 0.5 , X + Y = 0.7963643
cat("Cuantiles al 95%: X =", quantile_X_95, ", Y =", quantile_Y_95, ", X + Y =", quantile_Sum_XY_95_Negative, "\n")Cuantiles al 95%: X = 2.995732 , Y = 1.497866 , X + Y = 3.000276
cat("Valor esperado condicional de la cola al 95%: X =", 2, ", Y =", 2, ", X + Y =", mean_Sum_XY_Tail_95_Negative, "\n")Valor esperado condicional de la cola al 95%: X = 2 , Y = 2 , X + Y = 3.993384
# Imprimir los resultados
cat("Correlación Positiva:\n")Correlación Positiva:
cat("Varianzas: X =", var_X, ", Y =", var_Y, ", X + Y =", var_Sum_XY_Positivo, "\n")Varianzas: X = 1 , Y = 0.25 , X + Y = 2.139616
cat("Desviaciones Estándar: X =", sd_X, ", Y =", sd_Y, ", X + Y =", sd_Sum_XY_Positivo, "\n")Desviaciones Estándar: X = 1 , Y = 0.5 , X + Y = 1.462743
cat("Cuantiles al 95%: X =", quantile_X_95, ", Y =", quantile_Y_95, ", X + Y =", quantile_Sum_XY_95_Positivo, "\n")Cuantiles al 95%: X = 2.995732 , Y = 1.497866 , X + Y = 4.390258
cat("Valor esperado condicional de la cola al 95%: X =", 2, ", Y =", 2, ", X + Y =", mean_Sum_XY_Tail_95_Positivo, "\n")Valor esperado condicional de la cola al 95%: X = 2 , Y = 2 , X + Y = 5.877844
Cópula Marshall-Olkin
\[\hat{C}(u,v)=(1-u)(1-v)\textbf{mín}\{(1-u)^{-\alpha_{1}},(1-u)^{-\alpha_{1}}\}\]
Función de distribución conjunta
Primero notemos que esta sería la distribución conjunta de supervivencia de \(X\) y \(Y\):
\[S(x, y) = (e^{-\lambda x})(e^{-2\lambda y})\textbf{mín}\{(e^{-\lambda x})^{-\alpha_{1}},(e^{-2\lambda y})^{-\alpha_{1}}\}\]
Ahora mostraremos la distribución conjunta de \(X\) y \(Y\): \[H(x,y)= \]
Función de densidad conjunta
Arquimediana 1
\[C(u,v)=\phi^{-1}(\phi(u)+\phi(v))\]
donde \(\phi(u)=(-\ln u)^{\theta}\)
Arquimediana 2
\[C(u,v)=\phi^{-1}(\phi(u)+\phi(v))\]
donde \(\phi(u)=\ln(1-\theta\ln(u))\)
Arquimediana 3
\[C(u,v)=\phi^{-1}(\phi(u)+\phi(v))\]
donde \(\phi(u)=(1-\ln(u))^{\theta}-1\)
Caso Pareto
Tenemos las siguientes variables aleatorias:
\[X\sim \textbf{Pareto}(\alpha, 2\nu)\]
\[Y\sim \textbf{Pareto}(\alpha, \nu)\]
donde \(F_{x}(x)=1-\left(\frac{2\nu}{x+2\nu}\right)^{\alpha}\) y \(F_{y}(y)=1-\left(\frac{\nu}{y+\nu}\right)^{\alpha}\)
Cópula Gaussiana
\[C(u,v)=\Phi_{2}(\Phi^{-1}(u),\Phi^{-1}(v))\]
Cópula Marshall-Olkin
\[C(u,v)=(1-u)(1-v)\textbf{mín}\{(1-u)^{-\alpha_{1}},(1-u)^{-\alpha_{1}}\}\]
Arquimediana 1
\[C(u,v)=\phi^{-1}(\phi(u)+\phi(v))\]
donde \(\phi(u)=(-\ln u)^{\theta}\)
Arquimediana 2
\[C(u,v)=\phi^{-1}(\phi(u)+\phi(v))\]
donde \(\phi(u)=\ln(1-\theta\ln(u))\)
Arquimediana 3
\[C(u,v)=\phi^{-1}(\phi(u)+\phi(v))\]
donde \(\phi(u)=(1-\ln(u))^{\theta}-1\)